When PySAL was originally planned, the intention was to focus on the computational aspects of exploratory spatial data analysis and spatial econometric methods, while relying on existing GIS packages and visualization libraries for visualization of computations. Indeed, we have partnered with esri and QGIS towards this end.
However, over time we have received many requests for supporting basic geovisualization within PySAL so that the step of having to interoperate with an exertnal package can be avoided, thereby increasing the efficiency of the spatial analytical workflow.
In this notebook, we demonstrate several approaches towards geovisualization within a self-contained exploratory workflow. The idea here is the support quick generation of different views of your data to complement the statistical and econometric work in PySAL. Once your work has progressed to the publication stage, we point you to resources that can be used for publication quality output.
Contributors:
This document describes the main structure, components and usage of the mapping module in PySAL. The is organized around three main layers:
In [1]:
import numpy as np
import pysal as ps
import random as rdm
from pysal.contrib.viz import mapping as maps
%matplotlib inline
from pylab import *
This includes basic functionality to read spatial data from a file (currently only shapefiles supported) and produce rudimentary Matplotlib objects. The main methods are:
These methods all support an option to subset the observations to be plotted (very useful when missing values are present). They can also be overlaid and combined by using the setup_ax function. the resulting object is very basic but also very flexible so, for minds used to matplotlib this should be good news as it allows to modify pretty much any property and attribute.
In [2]:
shp_link = ps.examples.get_path('columbus.shp')
shp = ps.open(shp_link)
some = [bool(rdm.getrandbits(1)) for i in ps.open(shp_link)]
fig = figure()
base = maps.map_poly_shp(shp)
base.set_facecolor('none')
base.set_linewidth(0.75)
base.set_edgecolor('0.8')
some = maps.map_poly_shp(shp, which=some)
some.set_alpha(0.5)
some.set_linewidth(0.)
cents = np.array([poly.centroid for poly in ps.open(shp_link)])
pts = scatter(cents[:, 0], cents[:, 1])
pts.set_color('red')
ax = maps.setup_ax([base, some, pts], [shp.bbox, shp.bbox, shp.bbox])
fig.add_axes(ax)
show()
In [21]:
net_link = ps.examples.get_path('eberly_net.shp')
net = ps.open(net_link)
values = np.array(ps.open(net_link.replace('.shp', '.dbf')).by_col('TNODE'))
pts_link = ps.examples.get_path('eberly_net_pts_onnetwork.shp')
pts = ps.open(pts_link)
fig = figure()
netm = maps.map_line_shp(net)
netc = maps.base_choropleth_unique(netm, values)
ptsm = maps.map_point_shp(pts)
ptsm = maps.base_choropleth_classif(ptsm, values)
ptsm.set_alpha(0.5)
ptsm.set_linewidth(0.)
ax = maps.setup_ax([netc, ptsm], [net.bbox, net.bbox])
fig.add_axes(ax)
show()
In [22]:
maps.plot_poly_lines(ps.examples.get_path('columbus.shp'))
In [23]:
shp_link = ps.examples.get_path('columbus.shp')
values = np.array(ps.open(ps.examples.get_path('columbus.dbf')).by_col('HOVAL'))
types = ['classless', 'unique_values', 'quantiles', 'equal_interval', 'fisher_jenks']
for typ in types:
maps.plot_choropleth(shp_link, values, typ, title=typ)
Contributors:
In addition to using matplotlib, the viz module includes components that interface with the folium library which provides a Pythonic way to generate Leaflet maps.
In [24]:
import pysal as ps
import geojson as gj
from pysal.contrib.viz import folium_mapping as fm
First, we need to convert the data into a JSON format. JSON, short for "Javascript Serialized Object Notation," is a simple and effective way to represent objects in a digital environment. For geographic information, the GeoJSON standard defines how to represent geographic information in JSON format. Python programmers may be more comfortable thinking of JSON data as something akin to a standard Python dictionary.
In [25]:
filepath = ps.examples.get_path('south.shp')[:-4]
shp = ps.open(filepath + '.shp')
dbf = ps.open(filepath + '.dbf')
In [26]:
js = fm.build_features(shp, dbf)
Just to show, this constructs a dictionary with the following keys:
In [27]:
js.keys()
Out[27]:
In [28]:
js.type
Out[28]:
In [29]:
js.bbox
Out[29]:
In [30]:
js.features[0]
Out[30]:
Then, we write the json to a file:
In [31]:
with open('./example.json', 'w') as out:
gj.dump(js, out)
In [32]:
list(js.features[0].properties.keys())[:5]
Out[32]:
We can map these attributes by calling them as arguments to the choropleth mapping function:
In [15]:
import folium as fol
In [16]:
fm.choropleth_map('./example.json', 'FIPS', 'HR90')
Out[16]:
This produces a map using default classifications and color schemes and saves it to an html file. We set the function to have sane defaults. However, if the user wants to have more control, we have many options available.
There are arguments to change the classification scheme:
In [17]:
fm.choropleth_map('./example.json', 'FIPS', 'HR90', classification = 'Quantiles')
Out[17]: